// the compiled Rs-DeltaS relation
// created: March 15, 2021 by Yazhou SUN
// lasted modified: May 06, 2022

double Rs(double rs, bool opt);
double DS(double ds, bool opt);

struct tRs{
  string name;
  double rs, drs;
  double ds;
  tRs(const string &n, double rs_, double drsm, double drsp, double ds_, bool opt1 = true, bool opt2 = true){
    name = n;
    if(fabs(rs_-(drsm+drsp)/2.) > 0.3){
      cout << name << ": large discrepancy happened." << endl;
      cout << "rs_: " << rs_ << " drsm: " << drsm << " drsp: " << drsp << endl;
      getchar();
    }

    rs = Rs(rs_, opt1);
    drs = fabs(Rs(drsm, opt1)-Rs(drsp, opt1))/2.;
    ds = DS(ds_, opt2);
  } // end the constructor
  void Print(){
    cout << name << " rs: " << rs << "(" << drs << ")," << " ds: " << ds << endl;
  }
};

void rsdsCmpHF(){
  vector<tRs> vn, vp, ve;

  vn.push_back(tRs("Be10-n", 0.35, -0.7, 1.4, 11.1));
  vn.push_back(tRs("Ar46-n", 4.5, 1.4, 7.7, 12.55));
  vn.push_back(tRs("Li9-n", 4.75, 3.75, 5.85, 13.02));
  vn.push_back(tRs("Si38-n", 4.8, 4.0, 5.6, 9.27));
  vn.push_back(tRs("Si40-n", 5.0, 4.0, 6.1, 7.81));
  vn.push_back(tRs("Si36-n", 3.0, 2.45, 3.5, 10.8));
  vn.push_back(tRs("C15-n", 7.45, 6.4, 8.5, 6.51));
  vn.push_back(tRs("C19-n", 7.9, 5.8, 10.0, 2.45));
  vn.push_back(tRs("C20-n", 8.2, 6.4, 10.0, 2.80));
  vn.push_back(tRs("Ni57-n", -3.8, -4.3, -3.3, 11.75, 0,0));
  vn.push_back(tRs("C12-n", -3.55, -4.05, -3.0, 13.74, 0,0));
  vn.push_back(tRs("O16-n", -1.75, -2.5, -0.95, 10.75, 0,0));
  vn.push_back(tRs("Mg22-n", -7.2, -7.9, -6.4, 6.66, 0,0));
  vn.push_back(tRs("Ar34-n", -5.65, -7.95, -3.8, 6.62, 0,0));
  vn.push_back(tRs("Ca36-n", -10.0, -10.55, -9.4, 4.8, 0,0));
  vn.push_back(tRs("Ar32-n", -9.75, -10.55, -8.9, 3.16, 0,0));
  vn.push_back(tRs("S28-n", -8.2, -7.4, -8.9, 3.26, 0,0));
  vn.push_back(tRs("Si24-n", -6.15, -7.2, -5.1, 4.1, 0,0));
  vn.push_back(tRs("C10-n", -3.6, -4.05, -3.0, 4.4, 0,0));

  vp.push_back(tRs("O14-p", 2.2, 5.1, -0.7, 7.38));
  vp.push_back(tRs("Si24-p", 3.2, 2.2, 4.3, 7.86));
  vp.push_back(tRs("Ca36-p", 3.75, 2.7, 4.8, 8.5));
  vp.push_back(tRs("C9-p", 8.4, 7.1, 9.7, 10.8));
  vp.push_back(tRs("C9-p", 3.75, 2.2, 5.3, 11.05));
  vp.push_back(tRs("Be8-p", 4.8, 3.0, 6.65, 11.15));
  vp.push_back(tRs("S28-p", 6.4, 4.6, 8.2, 7.06));
  vp.push_back(tRs("C12-p", -3.7, -4.25, -3.2, 17.35, 1,0));
  vp.push_back(tRs("O16-p", 1.4, 0.4, 2.4, 15.25, 0,0));
  vp.push_back(tRs("Be10-p", -5.35, -5.9, -4.8, 7.37, 0,0));
  vp.push_back(tRs("Mg28-p", 0.6, 0.1, 1.1, 10.75, 0,0));
  vp.push_back(tRs("Si36-p", -6.15, -6.65, -5.6, 7.00, 0,0));
  vp.push_back(tRs("Si38-p", -8.45, -8.9, -7.9, 5.53, 0,0));
  vp.push_back(tRs("Si40-p", -8.45, -10.25, -6.65, 4.04, 0,0));

  ve.push_back(tRs("Ca40-e", -0.65, -1.95, 0.65, 14.76));
  ve.push_back(tRs("P31-e", 0.15, -1.1, 1.15, 18.95, 1,0));
  ve.push_back(tRs("O16-e", -1.1, -2.7, 0.9, 18.0, 1,0));
  ve.push_back(tRs("C12-e", -1.9, -0.9, -3.0, 17.5, 1,0));
  ve.push_back(tRs("Zr90-e", -3.0, -4.25, -1.7, 18.05, 1,0));
  ve.push_back(tRs("V51-e", -4.8, -5.8, -3.75, 17.7, 1,0));
  ve.push_back(tRs("Li7-e", 0.1, 1.6, -1.5, 13.9, 0,0));
  ve.push_back(tRs("Si30-e", -1.2, -2.5, 0.1, 13.8, 0,0));
  ve.push_back(tRs("Ca48-e", -2.25, -3.3, -1.2, 11.85, 0,0));
  ve.push_back(tRs("Pb208-e", -3.55, -4.8, -2.2, 15.35, 0,0));

  const int nn = vn.size(), np = vp.size(), ne = ve.size();
  double rs_n[nn], ers_n[nn], ds_n[nn], eds_n[nn];
  double rs_p[np], ers_p[np], ds_p[np], eds_p[np];
  double rs_e[ne], ers_e[ne], ds_e[ne], eds_e[ne];
  for(int i = 0; i < nn; i++){
    rs_n[i] = vn[i].rs; ers_n[i] = vn[i].drs; ds_n[i] = vn[i].ds; eds_n[i] = 0.;
    // vn[i].Print();
  }
  for(int i = 0; i < np; i++){
    rs_p[i] = vp[i].rs; ers_p[i] = vp[i].drs; ds_p[i] = vp[i].ds; eds_p[i] = 0.;
    // vp[i].Print();
  }
  for(int i = 0; i < ne; i++){
    rs_e[i] = ve[i].rs; ers_e[i] = ve[i].drs; ds_e[i] = ve[i].ds; eds_e[i] = 0.;
    // ve[i].Print();
  }
  TGraphErrors *gn = new TGraphErrors(nn, ds_n, rs_n, eds_n, ers_n);
  TGraphErrors *gp = new TGraphErrors(np, ds_p, rs_p, eds_p, ers_p);
  TGraphErrors *ge = new TGraphErrors(ne, ds_e, rs_e, eds_e, ers_e);
  TCanvas *c = new TCanvas("c", "Rs-DeltaS", 800, 600);
  gn->SetMinimum(0.15); gn->SetMaximum(1.3);
  gn->SetTitle("");
  TAxis *ax = gn->GetXaxis(), *ay = gn->GetYaxis();
  ax->SetRangeUser(-30, 24);
  ax->SetTitle("\\DeltaS (MeV)"); ay->SetTitle("R_{s}=\\sigma_{exp}/\\sigma_{th}");

  ax->SetLabelSize(0.05); ax->SetTitleSize(0.05);
  ay->SetLabelSize(0.05); ay->SetTitleSize(0.05);
  ax->SetTitleOffset(0.85); ay->SetTitleOffset(0.95);
  ax->CenterTitle(); ay->CenterTitle();

  double markersize = 2, linewidth = 1, col = 17;
  gn->SetMarkerStyle(20); gn->SetMarkerSize(markersize);
  gn->SetMarkerColor(17); gn->SetLineColor(col); gn->SetLineWidth(linewidth);
  gp->SetMarkerStyle(22); gp->SetMarkerSize(markersize);
  gp->SetMarkerColor(40); gp->SetLineColor(col); gp->SetLineWidth(linewidth);
  ge->SetMarkerStyle(21); ge->SetMarkerSize(markersize);
  ge->SetMarkerColor(45); ge->SetLineColor(col);  ge->SetLineWidth(linewidth);
  gStyle->SetEndErrorSize(9);
  TLegend *lg = new TLegend(0.105,0.105,0.45,0.3);
  lg->AddEntry(ge, "(e,e' p): \\DeltaS=S_{p}-S_{n}", "p");
  lg->AddEntry(gn, "n-removal: \\DeltaS=S_{n}-S_{p}", "p");
  lg->AddEntry(gp, "p-removal: \\DeltaS=S_{p}-S_{n}", "p");
  lg->SetBorderSize(0);

  // r0 = HF //
  TGraphErrors *gpHF = new TGraphErrors();
  gpHF->SetPoint(0, -15.358900, 0.872783); gpHF->SetPointError(0, 0., 0.135208); // O13
  gpHF->SetPoint(1, -18.552300, 0.700709); gpHF->SetPointError(1, 0., 0.054570); // O14
  gpHF->SetPoint(2, -2.411762, 0.634173); gpHF->SetPointError(2, 0., 0.058781); // O15
  gpHF->SetPoint(3, 0.681716, 0.428275); gpHF->SetPointError(3, 0., 0.024082); // O16
  gpHF->SetPoint(4, 9.813853, 0.504571); gpHF->SetPointError(4, 0., 0.089597); // O17
  gpHF->SetPoint(5, 10.955386, 0.288749); gpHF->SetPointError(5, 0., 0.021461); // O18
  gpHF->SetPoint(6, 13.597648, 0.431776); gpHF->SetPointError(6, 0., 0.069446); // O19
  gpHF->SetPoint(7, 12.215451, 0.506501); gpHF->SetPointError(7, 0., 0.026166); // O20
  gpHF->SetPoint(8, 16.812843, 0.423334); gpHF->SetPointError(8, 0., 0.028888); // O22
  gpHF->SetMarkerStyle(20); gpHF->SetMarkerSize(2);
  gpHF->SetMarkerColor(kMagenta); gpHF->SetLineColor(kMagenta);
  gpHF->SetLineWidth(linewidth+2);
  // p2p results //
  TGraphErrors *g2pHF = new TGraphErrors();
  g2pHF->SetPoint(0, -15.358900, 0.689954); g2pHF->SetPointError(0, 0., 0.152792); // O13
  g2pHF->SetPoint(1, -18.552300, 0.523866); g2pHF->SetPointError(1, 0., 0.074253); // O14
  g2pHF->SetPoint(2, -2.382810, 0.488594); g2pHF->SetPointError(2, 0., 0.077989); // O15
  g2pHF->SetPoint(3, 0.748679, 0.420615); g2pHF->SetPointError(3, 0., 0.040745); // O16
  g2pHF->SetPoint(4, 9.813900, 0.424770); g2pHF->SetPointError(4, 0., 0.040864); // O17
  g2pHF->SetPoint(5, 11.119888, 0.383515); g2pHF->SetPointError(5, 0., 0.046754); // O18
  g2pHF->SetPoint(6, 17.770223, 0.344451); g2pHF->SetPointError(6, 0., 0.036975); // O21
  g2pHF->SetPoint(7, 16.811993, 0.299550); g2pHF->SetPointError(7, 0., 0.039906); // O22
  g2pHF->SetPoint(8, 20.973753, 0.298327); g2pHF->SetPointError(8, 0., 0.077410); // O23
  g2pHF->SetMarkerStyle(22); g2pHF->SetMarkerSize(2);
  g2pHF->SetMarkerColor(kBlue); g2pHF->SetLineColor(kBlue);
  g2pHF->SetLineWidth(linewidth+2);
  // p2p-Rs by GSI //
  TGraphErrors *g2pgsi = new TGraphErrors();
  g2pgsi->SetPoint(0, -18.552300, 0.68); g2pgsi->SetPointError(0, 0., 0.07); // O14
  g2pgsi->SetPoint(1, 0.748679, 0.70); g2pgsi->SetPointError(1, 0., 0.05); // O16
  g2pgsi->SetPoint(2, 9.813900, 0.65); g2pgsi->SetPointError(2, 0., 0.05); // O17
  g2pgsi->SetPoint(3, 17.770223, 0.58); g2pgsi->SetPointError(3, 0., 0.04); // O21
  g2pgsi->SetPoint(4, 20.973753, 0.62); g2pgsi->SetPointError(4, 0., 0.13); // O23
  g2pgsi->SetMarkerStyle(23); g2pgsi->SetMarkerSize(2);
  g2pgsi->SetMarkerColor(kGreen); g2pgsi->SetLineColor(kGreen);
  g2pgsi->SetLineWidth(linewidth+2);


  // TGraphErrors *gnHF = new TGraphErrors();
  // gnHF->SetPoint(0, -17.903, 0.854988); gnHF->SetPointError(0, 0., 0.086473); // C16
  // gnHF->SetPoint(1, -17.591669, 0.865110); gnHF->SetPointError(1, 0., 0.091392); // C15
  // gnHF->SetPoint(2, -10.610312, 0.665784); gnHF->SetPointError(2, 0., 0.061324); // C14
  // gnHF->SetPoint(3, 18.3031, 0.278293); gnHF->SetPointError(3, 0., 0.030406); // C16-p
  // gnHF->SetMarkerStyle(20); gnHF->SetMarkerSize(1);
  // gnHF->SetMarkerColor(kYellow); gnHF->SetLineColor(kMagenta);
  // gnHF->SetLineWidth(linewidth+2);
  // // r0 = 1.25fm //
  // TGraphErrors *gnNew = new TGraphErrors();
  // gnNew->SetPoint(0, -17.879, 1.057552); gnNew->SetPointError(0, 0., 0.106960); // C16
  // gnNew->SetPoint(1, -17.586, 0.760564); gnNew->SetPointError(1, 0., 0.080348); // C15
  // gnNew->SetPoint(2, -10.7592, 0.726234); gnNew->SetPointError(2, 0., 0.066892); // C14
  // gnNew->SetPoint(3, 18.3031, 0.32339059); gnNew->SetPointError(3, 0., 0.029724); // C16-p
  // gnNew->SetMarkerStyle(20); gnNew->SetMarkerSize(markersize+0.2);
  // gnNew->SetMarkerColor(kRed); gnNew->SetLineColor(kRed);
  // gnNew->SetLineWidth(linewidth+2);
  // // r0 = 1.15fm //
  // const double dSOfst = 1.; // ds-offset
  // TGraphErrors *gnNew1 = new TGraphErrors();
  // gnNew1->SetPoint(0, -17.885510-dSOfst, 1.133221); gnNew1->SetPointError(0, 0., 0.114613); // C16
  // gnNew1->SetPoint(1, -17.650128-dSOfst, 0.796219); gnNew1->SetPointError(1, 0., 0.084115); // C15
  // gnNew1->SetPoint(2, -10.754304-dSOfst, 0.782902); gnNew1->SetPointError(2, 0., 0.072112); // C14
  // gnNew1->SetPoint(3, 18.3031-dSOfst, 0.35245669); gnNew1->SetPointError(3, 0., 0.0323951); // C16-p
  // gnNew1->SetMarkerStyle(23); gnNew1->SetMarkerSize(2);
  // gnNew1->SetMarkerColor(1); gnNew1->SetLineColor(1);
  // gnNew1->SetLineStyle(7); gnNew1->SetLineWidth(linewidth+1);
  // // r0 = 1.35fm //
  // TGraphErrors *gnNew2 = new TGraphErrors();
  // gnNew2->SetPoint(0, -17.872253+dSOfst, 0.987041); gnNew2->SetPointError(0, 0., 0.099829); // C16
  // gnNew2->SetPoint(1, -17.520026+dSOfst, 0.726037); gnNew2->SetPointError(1, 0., 0.076700); // C15
  // gnNew2->SetPoint(2, -10.762913+dSOfst, 0.673872); gnNew2->SetPointError(2, 0., 0.062069); // C14
  // gnNew2->SetPoint(3, 18.3031+dSOfst, 0.29669815); gnNew2->SetPointError(3, 0., 0.0272702); // C16-p
  // gnNew2->SetMarkerStyle(22); gnNew2->SetMarkerSize(2);
  // gnNew2->SetMarkerColor(4); gnNew2->SetLineColor(4);
  // gnNew2->SetLineStyle(2); gnNew2->SetLineWidth(linewidth+1);

  TLegend *lgNew = new TLegend(0.505,0.65,0.895,0.86);
  // lgNew->AddEntry(gnNew1, "r_{0}=1.15 fm", "pl");
  // lgNew->AddEntry(gnNew, "r_{0}=1.25 fm", "pl");
  // lgNew->AddEntry(gnNew2, "r_{0}=1.35 fm", "pl");
  // lgNew->AddEntry(gnHF, "r_{0} from HF", "pl");
  lgNew->AddEntry(gpHF, "r_{0} from HF_O-p@C-TA", "pl");
  lgNew->AddEntry(g2pHF, "r_{0} from HF_O-(p,2p)", "pl");
  lgNew->AddEntry(g2pgsi, "r_{0} from HF_O-(p,2p)-GSI", "pl");
  lgNew->SetBorderSize(0);

  const int n = 2;
  double xh[n] = {-30., 30.}, yh[n] = {1.23, 0.4};
  TGraph *gHatch = new TGraph(2, xh, yh);
  gHatch->SetLineColor(kWhite);
  gHatch->SetLineWidth(-3200);
  gHatch->SetFillStyle(3004);
  gHatch->SetFillColor(1);

  gn->Draw("Ap");
  ge->Draw("samep");
  gp->Draw("samep");
  gHatch->Draw("sameC");
  // gnNew->Draw("samepl");
  // gnNew1->Draw("samepl");
  // gnNew2->Draw("samepl");
  // gnHF->Draw("samepl");
  TF1 *fpol1p = new TF1("fpol1p", "pol1", -30., 25.);
  TF1 *fpol12p = new TF1("fpol12p", "pol1", -30., 25.);
  TF1 *fpol12pgsi = new TF1("fpol12pgsi", "pol1", -30., 25.);
  fpol1p->SetLineColor(gpHF->GetLineColor());
  fpol12p->SetLineColor(g2pHF->GetLineColor());
  fpol12pgsi->SetLineColor(g2pgsi->GetLineColor());
  gpHF->Fit("fpol1p");
  g2pHF->Fit("fpol12p");
  g2pgsi->Fit("fpol12pgsi");
  gpHF->Draw("samep");
  g2pHF->Draw("samep");
  g2pgsi->Draw("samep");
  lg->Draw();
  lgNew->Draw();

  // TPaveText *pt1 = new TPaveText(-24.72224,0.8263144,-20.38554,0.9418366,"br");
  // pt1->SetBorderSize(0); pt1->SetFillColor(0); pt1->SetFillColorAlpha(0, 1);
  // pt1->SetTextFont(42); pt1->AddText("^{16}C-n"); pt1->SetTextColor(2); pt1->Draw();
  // TPaveText *pt2 = new TPaveText(-23.94021,0.5494009,-19.60351,0.6649231,"br");
  // pt2->SetBorderSize(0); pt2->SetFillColor(0); pt2->SetFillColorAlpha(0, 1);
  // pt2->SetTextFont(42); pt2->AddText("^{15}C-n"); pt2->SetTextColor(2); pt2->Draw();
  // TPaveText *pt3 = new TPaveText(-16.11993,0.4066971,-11.78322,0.5222192,"br");
  // pt3->SetBorderSize(0); pt3->SetFillColor(0); pt3->SetFillColorAlpha(0, 1);
  // pt3->SetTextFont(42); pt3->AddText("^{14}C-n"); pt3->SetTextColor(2); pt3->Draw();
  // TPaveText *pt4 = new TPaveText(19.35902,0.367296,23.69813,0.4827474,"br");
  // pt4->SetBorderSize(0); pt4->SetFillColor(0); pt4->SetFillColorAlpha(0, 1);
  // pt4->SetTextFont(42); pt4->AddText("^{16}C-p"); pt4->SetTextColor(2); pt4->Draw();
  // for(int i = 13; i <= 23; i++){
  //   double j = i-13;
  //   double x0 = -24.7222+j*4., x1 = -20.3855+j*4.;
  //   double y0 = 0.8263-j*0.08, y1 = 0.9418-j*0.08;
  //   TPaveText *pt = new TPaveText(x0, y0, x1, y1, "br");
  //   pt->SetBorderSize(0); pt->SetFillColor(0); pt->SetFillColorAlpha(0, 1);
  //   pt->SetTextFont(42); pt->SetTextColor(2);
  //   char txt[128]; sprintf(txt, "^{%d}O", i); pt->AddText(txt); pt->Draw();
  // } // end for over i

  Int_t ci;      // for color index setting
  TPaveText *pt = new TPaveText(-14.79058,0.9288457,-12.41888,0.9909163,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  TText *pt_LaTex = pt->AddText("^{13}O");
  pt->Draw();

  pt = new TPaveText(-21.58708,0.4252269,-19.21539,0.4872976,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{14}O");
  pt->Draw();

  pt = new TPaveText(-3.569256,0.7059556,-1.197558,0.7680263,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{15}O");
  pt->Draw();

  pt = new TPaveText(-1.055964,0.2926214,1.315734,0.3546921,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{16}O");
  pt->Draw();

  pt = new TPaveText(7.067987,0.5296185,9.439686,0.5916892,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{17}O");
  pt->Draw();

  pt = new TPaveText(8.501626,0.196694,10.87332,0.2587647,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{18}O");
  pt->Draw();

  pt = new TPaveText(12.46626,0.2926214,14.83795,0.3546921,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{19}O");
  pt->Draw();

  pt = new TPaveText(10.83793,0.5437255,13.20962,0.6057961,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{20}O");
  pt->Draw();

  pt = new TPaveText(17.52824,0.2333722,19.89994,0.2954428,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{21}O");
  pt->Draw();

  pt = new TPaveText(15.36893,0.461905,17.74063,0.5239757,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{22}O");
  pt->Draw();

  pt = new TPaveText(20.07693,0.3885488,22.44863,0.4506195,"br");
  pt->SetBorderSize(0);
  pt->SetFillColor(0);
  pt->SetTextColor(2);
  pt->SetTextFont(42);
  pt_LaTex = pt->AddText("^{23}O");
  pt->Draw();


  c->Print("rsCmpHF.pdf");
  // c->Print("rsCmp.eps");
} // end the main function



// -----------
double Rs1(double rs){
  static const double kRs = 0.8/(11.+9.82), bRs = 0.76-2.2*kRs;
  return kRs*rs+bRs;
}
double Rs2(double rs){
  static const double kRs = 0.8/(11.+9.82), bRs = 0.5+3.85*kRs;
  return kRs*rs+bRs;
}
double DS1(double ds){
  static const double kds = 16./10.5, bds = -10.-13.02*kds;
  return kds*ds+bds;
}
double DS2(double ds){
  static const double kds = -16./10.5, bds = -3.-17.7*kds;
  return kds*ds+bds;
}
double Rs(double rs, bool opt){
  if(opt) return Rs1(rs);
  return Rs2(rs);
}
double DS(double ds, bool opt){
  if(opt) return DS1(ds);
  return DS2(ds);
}
